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.*-< \ Taylor-Couette flow with independently rotating inner (i) & outer (o) cylinders is ex- 

plored numerically and experimentally to determine the effects of the radius ratio r\ on the 
system response. Numerical simulations reach Reynolds numbers of up to Re% = 9.5 • 10 3 
O ,■ and Re = 5-10 3 , corresponding to Taylor numbers of up to Ta = 10 8 for four dif- 

ferent radius ratios r\ = ri/r between 0.5 and 0.909. The experiments, performed in 
the Twente Turbulent Taylor-Couette (T 3 C) setup, reach Reynolds numbers of up to 
Re, = 2 • 10 6 and Re = 1.5 • 10 6 , corresponding to Ta = 5 ■ 10 12 for r) = 0.714 - 0.909. 
Effective scaling laws for the torque J u (Ta) are found, which for sufficiently large driving 
Ta arc independent of the radius ratio rj. As previously reported for r\ = 0.714, optimum 
transport at a non-zero Rossby number Ro = rj|wj — uj \/[2(r — r.;)oj ] is found in both 
experiments and numerics. Ro op t is found to depend on the radius ratio and the driving 
of the system. At a driving of about Ta ~ 10 10 , Ro opt saturates to an asymptotic t\- 
dependent value. Theoretical predictions for the asymptotic value of Ro op t are compared 
to the experimental results, and found to differ notably. Furthermore, the local angular 
velocity profiles from experiments and numerics are compared, and a link between a flat 
bulk profile and optimum transport for all radius ratios is reported. 
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1. Introduction 

Taylor-Couette (TC) flow consists of the flow between two coaxial cylinders which are 
independently rotating. A schematic drawing of the system can be seen in Fig. [1] The 
rotation difference between the cylinder shears the flow thus driving the system. This 
rotation difference has been traditionally expressed by two Reynolds numbers, the inner 
cylinder Re h = riUJid/v, and the outer cylinder Re = r o L0 o d/v Reynolds numbers, where 
n and r are the radii of the inner and outer cylinder, respectively, w, and u> the inner and 
outer cylinder angular velocity, d = r a — r% the gap width, and v the kinematic viscosity. 
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Figure 1: Schematic of the Taylor- Coucttc system. The system consists of two coaxial 
cylinders, which have an inner cylinder radius of n and an outer cylinder radius of r . 
Both cylinders are of length L. The inner cylinder rotates with an angular velocity Wi 
and the outer cylinder rotates with an angular velocity of ui . 

The geometry of TC is characterized by two nondimensional parameters, namely the 
radius ratio 77 = Tj/r and the aspect ratio T = L/d. 

Instead of taking Rei and Re , the driving in TC can alternatively be characterized 
by the Taylor Ta and the rotation rate, also called the Rossby Ro number. The Taylor 
number can be seen as the non-dimensional forcing (the differential rotation) of the 
system defined as Ta = a(r — ri) 2 (r + r.;) 2 (cj — uji) 2 /(4i/ 2 ), or 



Here 



r t/ r t With r a = ( r o 



Ta=(r h a d 2 /ririv 2 )(u ~u Ji y. 
■Ti)/2 the arithmetic and r g = 



(1.1) 
lr Ti the geometric mean 



radii. The Rossby number is defined as: 



Ro 



2uj d 



(1.2) 



and can be seen as a measure of the rotation of the system as a whole. Ro < corresponds 
to counterrotating cylinders, and Ro > to corotating cylinders. 

TC is among the most investigated systems in fluid mechanics, mainly owing to its 
simplicity as an experimental model for shear flows. TC is in addition a closed system, so 
global balances which relate the angular velocity transport to the energy dissipation can 
be obtained. Specifically, in Eckhardt, Grossmann & Lohse (2007) (from now on referred 
to as EGL 2007), an exact relationship between the global parameters and the volume 
averaged energy dissipation rate was derived. This relationship has an analogous form to 
the one which can be obtained for Rayleigh-Bcnard (RB) flow, i.e. a flow in which heat 
is transported from a hot bottom plate to a cold top plate. 

TC and RB flow have been extensively used to explore new concepts in fluid me- 



Optimal Taylor- Couette flow: Radius ratio dependence 



chanics. Instabilities (IPfister fe Rehberd [l98lt IPfister et g]]|l988L IChandrasekhail 1 1 98 ll 
Drazin fe Reidl 1981t Bussdll967 ), nonlinear dynamic s and chaos (ILorena Il963t lAhlersI 
1974 iBehringerl Il985[ iDominguez-Lerma et al. 1986 ; Strogata 1994 ) patte rn forma- 



tion (lAndereck et aZj|1986t Cross fe Hohenbcrg||l993l: iBodenschatz et al\ 200(J). and tur 



bulence (ISiggial 11994; iGrossmann fe Lohsd 12000: lKadanofll2001t lLathrop et dj 119926 : 
Ahlers et all l2009t lLohse fe Xiall2010l) have been studied in both TC and RB and both 
numerically and experimentally. The main reasons behind the popularity of these sys- 
tems are, in addition to the fact that they are closed systems, as mentioned previously, 
their simplicity due to the high amount of symmetries present. It is also worth noting 
that plane Couette flow is the limiting case of TC when the radius ratio r\ = 1. 

Experi mental invest igati ons of TC ha ve a long history, dating back to the seminal 
works by Wendt J 1933 ) and Tavlorl ( 19361). Experimental work continued during the yea rs 



(ISmith fe Townsendlll98"3 lAndereck et al.lll986l : lTong et aZ.lll99d;lLathrop et ali\l992^!ia]: 
Lewis fe Swinnevlll999l ; Ivan Gils et al\\2011a\ \a; IPaoletti fe Lathrop|l201ll : IHuisman et al 



201261 ) at low and high Ta numbers for different ratios of the rotation frequencies 
a = —uj /uji. a is positive for counter-rotation and negative for co-rotation, —a = /j, 
another measure used for the ratio of rotation frequencies. This work has been comple- 
ment ed by numerical simulations, not only in the regime of pure inner cylinder rota- 
tion (JFasel fc Boozlll984 ICoughlin fc Marcuslll996t |Pong||2007l 120081: IpTrro fc Quadrfol 
20081 ). but also for eige nvalue study ( Gebhardt fc GrossmanrJ 19931) . and counter-rotation 
at fixed a (|Dongll2008l ). Recently (JBrauckmann fc Eckhardll2013ol : lOstilla et aLll20136 ; 



Heise et al 



20131 ). simulations have also explored the effect of the outer cylinder rotation 



on the system. ^_^ 

Th e recent experime nts ( van Gils et a0l2011a A Paoletti fc Lathropl2011 : lMerbold et al 
20131 ) and simulations ( Brauckmann fc Eckhardtll2Q13al : l0stilla et aZ.ll20136i ) have shown 
that at fixed Ta an optimal angular momentum transport is obtained at non-zero a op t, 



and that the location of this maximum a opt varies with Ta. However, both experi- 
ments and simulations have been restricted to two radius ratios, namely r\ = 0.5 and 
?/ = 0.714. The same radius ratios were also used for studies carried out on scaling 
laws of the torque and the "wind" of turbulence at highly turbulent Taylor numbers 
dLewis fc Swinnev 1999 : Paoletti fc Lathrop||2011tlvan Gils et QJ.Il20116t IHuisman et al. 



201261 : iMerbold et q/.ll2013l ). Up to now, it is not clear how the radius ratio affects the 



scaling laws of the system response or the recently found phenomena of optimal transport 
as a function of Ta. 

Two suggestions were made to account for the radius ratio dependence of optimal trans- 
port. Van Gils et. al (2011b) wondered whether the optimal transport in general lies in 
or at least close to the Voronoi boundary (meaning a line of equal distance) of the Esser- 
Grossmann stability lines in the (Re , Rei) phase space as it does for r\ = 0.714. However , 
this bisector value does not give the co rrect optim al transport for r\ = 0.5 (IMerbold et al 



2013tlBrauckmann fc EckhardtH20i36h . Therefore ISrauckmann fc Eckhardd (|2013al ) pro 



posed a dynamic extension of the Esser-Grossmann instability theory. This model cor 
rectly gives the observed optimal transport (withi n experimental error bars) between n = 
0.5 and n = 0.714 for three experimental data sets (IWenddll933l:IPaoletti fc Lathrop|2011 



rtWen< 
; (|Brai 



but 



Ivan Gils et aLll20116l ) and one numerical data set (JBrauckmann fc Eckhardtll2013? 
it is not clear how it performs outside the 77-range [0.5, 0.714]. 

In this paper, we study the following questions: how does the radius ratio r\ affect the 
flow? How are the scaling laws of the angular momentum transport affected? What is the 
role of the geometric parameter called pseudo-Prandtl number a introduced in EGL2007? 
Can the effect of the radius ratio be interpreted as a kind of non-Oberbcck-Boussinesq 
effect, analogous to this effect in Raylcigh-Bcnard flow? Finally, arc the predictions and 
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insights of Ivan Gils et all ( 2011&I ) , lOstilla et all ( 2013M ) and iBrauckmann fe Eckhardtl 
(|20136h on the optimal transport also valid for other values of 77? 

In order to answer these questions, both direct numerical simulations (DNS) and 
experiments have been undertaken. Numerical simulations, with periodic axial bound- 
ary conditions, have been performed using the finite-difference code previously used in 
Ostilla et all (J2013M) . In these simulations, three more values of 77 have been investigated: 
one in which the gap is larger (77 = 0.5), and two in which the gap is s maller (77 = 0.833 
and 0.909). With the previous simulations from lOstilla et al\ ( 2013a) at r\ = 0.714, a 
total of four radius ratios has been analyzed. 

Figure [2] shows the (Ta,l/Ro) parameter space explored in the simulations for the 
four selected values of the radius ratio 77. A higher density of points has been used in 
places where the global response (Nu u , Re w ) of the flow shows more variation with the 
control parameters Ta and l/Ro. A fixed aspect ratio of T = 2ir has been taken for all 
simulations, and axially periodic boundary conditions have been used. These sim ulations 
have the same upper bounds of Ta (or Rei) as the ones of lOstilla et al\ (|20136l ). 

In addition to these simulations, experiments have been performed with the Twcntc 
Turbulent Taylor-Couette (T 3 C) facility, with which we achieve larger Ta numbers. De- 
tails of the setup are given in Ivan Gils et all ( 2011qf ). Once again, four values of 77 have 
been investigated, but, due to experimental constraints, we have been limited to investi- 
gate only smaller gap widths, i.e. values 77 ^ 0.714, namely 77 = 0.714, 0.769, 0.833 and 
0.909. The experimentally explored parameter space are shown in Fig. [3] 

The structure of the paper is as follows. In sections 2 and 3, we start by describing 
the numerical code and the experimental setup, respectively. In section 4, the global 
response of the system, quantified by the non-dimensionalized angular velocity current 
Nu^ is analyzed. To understand the global response, we analyze the local data which 
can be obtained from the DNS simulations in section 5. Angular velocity profiles in the 
bulk and in the boundary layers are analyzed and related to the global angular velocity 
optimal transport. We finish in section 6 with a discussion of the results and an outlook 
for further investigations. 



2. Numerical method 

In this section, the used numerical method is explained in detail. The rotating frame in 
which the Navier-Stokes equations are solved and the employed non-dimensionalizations 
are introduced in the first section. This is followed by detailing the spatial resolution 
checks which have been performed. 

2.1. Code description 

The employed code is a finite difference code, which solves the Navier-Stokes equations in 
cylindrical coordinates. A second-order spatial discretization is used, and the equations 
are advanced in time by a fractional time integration method. This code is based on the 
so-cal led Verzicco-code, whose numerical algorithms are detailed in IVerzicco fc Orlandi 
( 19961 ). A combination of MPI and OpenMP directives are used to achieve large scale 
parallelization. This code has been extensivel y used for Rayleigh-Bcnard flow; for re- 
cent sim ulations see lStevens et all ( 2010ll201ll ). In the context of TC flow, lOstilla et al. 
( 201361 ) already validated the code for 77 = 0.714. 



The flow was simulated in a rotating frame, which was chosen to rotate with £7 = 0J o e z . 
This was done in order to simplify the boundary conditions. In that frame the outer 
cylinder is stationary for any value of a, while the inner cylinder has an azimuthal velocity 
of ug(r = n) = Ti{bj\ — u'), where the I superscript denotes variables in the lab frame, 
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Figure 2: Control parameter phase space which was numerically explored in this paper 
in the (Ta, 1/Ro) representation. From top-left to bottom-right: 77 = 0.5, 0.714, 0.833 
and 0.909. V = 2tt was fixed, and axial periodicity was employed. The grey-shaded 
area signals boundary conditions for which the angular momentum L — r 2 ui of the outer 
cylinder (L Q ) is larger than the angular momentum of the inner cylinder (Li). This causes 
the flow to have an overall transport of angular momentum towards the inner cylinder. 
In this region, the Rayleigh stability criterium applies, which states that if dL/dr > 
the flow is linearly stable to axisymmetric perturbations. 



while no superscript denotes variables in the rotating frame. We then choose the inner 
cylinder rotation rate in the rotating frame as the characteristic velocity of the system 
U = \ug(ri)\ = ri\uji —lu \ and the characteristic length scale d to non-dimensionalize the 
equations and boundary conditions. 

It is worth noting that the characteristic velocity U can be written as 



U = (v/d) ■ [8t7 2 /(1 + rjf] 1 / 2 ■ To 1 ' 2 . 



(2.1) 



This means that up to a geometric factor depending on rj, and which for our 77's is between 
0.6 and 1, the characteristic velocity Reynolds number Re = Udjv is simply Ta 1 / 2 . 

Using this non-dimensionalization, the inner cylinder velocity boundary condition sim- 
plifies to: iig(r = n) = sgn(uji — u) ). In this paper, LOi — lo is always positive. Thus, in 
this rotating frame the flow geometry is simplified to a pure inner cylinder rotation with 
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Figure 3: Control parameter phase space which was explored in experiments in the {Ta, 
Ro^ 1 ) representation for r\ = 0.716 (top-left), rj = 0.769 (top-right), r) = 0.833 (bottom- 
left) and 77 = 0.909 (bottom-right). 



the boundary condition ue( r «) = 1- The outer cylinder's effect on the flow is felt as a 
Coriolis force in this rotating frame of reference. The Navicr-Stokcs equations then read: 



Kl/ > V 2 u + Ro~ 1 e z x it 



du „ - „ - „ 

— - + u ■ \/u = — Vp 
dt \Ta 

where Ro was defined previously in Eq. 11.21 and f(r]) is 

(1 + r/) 3 



m 



8rf 



(2.2) 



(2-3) 



It is useful to continue the non-dimensionalization by defining the normalized radius 
f = (r — Ti)j d and the normalized height z = z/d. We define the time- and azimuthally- 
averaged velocity field as: 

u(r,z) = {u(0,r,z,t))g, t , (2.4) 

where (<p(xi, x%, ..., x n )) Xi indicates averaging of the field 4> with respect to Xi- 

To quantify the torque in the system, we first note that the angular velocity current 

r = r 3 ((uV) e , z>t - vd r {J) e , z ,t) (2.5) 

is conserved, i.e. independent on the radius r (EGL 2007). J u represents the current of 
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angular velocity from the inner cylinder to the outer cylinder (or vice versa). The first 
term is the convective contribution to the transport, while the second term is the diffusive 
contribution. 

In the state with the lowest driving, and ignoring end plate effects, a laminar, time 
independent velocity field which is purely azimuthal, u g (r) = Ar+B/r, with u r = u z = 0, 
is induced by the rotating cylinders. This laminar flow produces an angular velocity 
current J^, which can be used to nondimcnsionalize the angular velocity current and a 
non-zero dissipation rate e n? o. 

Nu u = ^j. (2.6) 

Nu^ can be seen as an angular velocity "Nusselt" number. 

When J" , and therefore Nu u are calculated numerically, the values will depend on the 
radial position, due to finite time averaging. We can define Aj to quantify this radial 
dependence as: 

max(J w (r)) — min( J u (r)) max.(Nu u (r)) — min(Nu u (r)) 
J = WWr = (Nu u (r)) r ( ' 7) 

which analytically equals zero but will deviate when calculated numerically. 

The convective dissipation per unit mass can be calculated either from its definition 
as a volume average of the local dissipation rate for an incompressible fluid, 

eu = ei = ^((dju j + d i jUi f) V!t , (2.8) 

or a global balance can also be used. The exact relationship (EGL 2007) 

4 - <o - ^<r- 2 Ta(N Uul - 1) , (2.9) 

where e u ,o is the volume averaged dissipation rate in the purely azimuthal laminar flow, 
links the volume averaged dissipation to the global driving Ta and response Nu^. 

This link can be and has been used for code validation and for checking spatial reso- 
lution adequatencss. The volume averaged dissipation can be calculated from both (|2.8|) 
and (|2.9j) and later checked for sufficient agreement. We define the quantity A £ as the rel- 
ative difference between the two ways of numerically calculating the dissipation, namely 
either via Nu u with cq. 12.91 or directly from the velocity gradients, eq. 12.81 

v 3 d- 4 <j- 2 Ta(N UuJ - 1) + e„, - § ((fifed + d l .u t f) 
A - - - ' — - - -lii (2 10) 

§<(#* + #*)% ' 

A e is equal to analytically, but will deviate when calculated numerically. The deviation 
of A, / and A c from zero is an indication of the adequateness of the resolution. 

We would like to emphasize that the requirement for Aj < 0.01 is much stricter than 
torque balan c e, whic h can simply be expressed as Nu u (n) = Nu u (r )- As analyzed in 



Ostilla et al\ (J20136I ). a value of less than 1% for A,/ and about 1% for A c is linked 
to grid adequateness at the Taylor number simulated. To ensure convergence in time, 
the time-averages of the Nusselt number and the energy dissipation calculated locally 
(equation 12. 8p were also checked to converge in time within 1%. 
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2.2. Resolution checks 

Spatial resolution checks were performed in two ways. First, as mentioned previously, the 
values of Aj and A e were checked. As an additional check, simulations at selected values 
of Ta were performed at a higher resolution. As the explored parameter space is large, 
these checks were performed only for the highest value of Ta simulated for the grid size. 
A lower driving of the flow for the same grid size is expected to have a smaller error due 
to spatial discretization, as spatial discretization errors increase with increased Re, and 
thus increased Ta. 

Concerning the temporal resolution there are numerical and physical constraints; the 
former requires a time step small enough to keep the integration scheme stable and 
this is achieved by using an adaptive time step based on a Courant-Frcdcrich-Lewy 
(CFL) criterium. The 3 rd -order Rungc Kutta time-marching algorithm allows for a 
CFL of up to y/3, but this can be reduced due to the implicit factorization of the 
viscous terms. For safety, the maximum CFL has been taken as 1.4. From the physical 
point of view, the time step size must also be small enough to properly describe the 
fast dynamics of the smallest flow scale which is the Kolmogorov scale. Although the 
time step size should be determined by the most restrictive among the two criteria 
above, our experience suggests that as long as the CFL number criterion is satisfied, 
which guarantees numerical stability, the results become insensitive to the time step size 
and all the flow scales are adequately described temporally. Direct confirmation of this 
statement can be found in lOstilla et al\ ( 2013a ). 



The results for r\ = 0.5, 0.833, and 0.909 are presented in Table [T] Uniform discretiza- 
tion was used in azimuthal and axial directions. In the radial direction, points were 
clustered near the walls by using hyperbolic tangent-type clustering, or a clipped Cheby- 
chev type clustering for higher values of Ta. A table including the results for the spatial 
resolution tests at 77 = 0.714 can be found in lOstilla et al\ ( 2013a ). 



3. Experimental setup 

The Twentc Turbulent Taylor-Couette (T 3 C) apparatus has been built to obtain high 



Ta nu mbers. It has been described in detail in Ivan Gils et all ( 2011af ) and Ivan Gils et 



(|2012l ). The inner cylinder with outside radius r, = 0.200 m consists of three sections. 
The total height of those axially stacked sections is L = 0.927 m. We measure the torque 
only on the middle section of the inner cylinder, which has a height of L m id = 0.536 m, to 
reduce end effects in our torque measurements. The transparent outer cylinder is made 
of acrylic and has an inside radius of r = 0.279 m. We vary the radius ratio by reducing 
the diameter of the outer cylinder by adding a 'filler' that is fixed to the outer cylinder 
and sits between the inner and the outer cylinder, effectively reducing r while keeping 
Ti fixed. We have 3 fillers giving us 4 possible outer radii: r = 0.279 m (without any 
filler), 0.26 m, 0.24 m, and 0.22 m, giving experimental access to n = 0.716, 0.769, 0.833, 
and 0.909, respectively. Note that by reducing the outer radius, we not only change 77, 
but also change T = L/(r - n) from T(r) = 0.716) = 11.68 to r(r? = 0.909) = 46.35. 

For high Ta the heating up of the system becomes apparent and it has to be actively 
cooled in order to keep the temperature constant. We cool the working fluid (water) 
from the top and bottom end plates and maintain the temperature constant during any 
experiment. The setup has been constructed in such a way that we can rotate both 
cylinders independently while keeping the setup cooled. 

As said before, we measure the torque on the middle inner cylinder. We do this by 
measuring the torque that is transferred from the axis to the cylinder by using a load-cell 
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Table 1: Resolution tests for T = 2ir and r; = 0.5, 0.833 and 0.909. The first column 
displays the radius ratio, the second column displays the Taylor number, the third column 
displays the resolution employed, the fourth column the calculated Nu u , the fifth column 
and sixth colums the relative discrepancies Aj and A e , and the last column the 'case': 
( R)esolved and ( E)rror reference. A e is positive, and exceeds the 1% threshold reported 
in lOstilla et all ( 2013M ) for some cases at the largest rj, but even so resolution appears 
to be sufficient as variations of Nu^ are small. 



that is inside the aforementioned cylinder. Torque measurements are performed using 
a fixed procedure: the inner cylinder is spun up to its maximum rotational frequency 
of 20 Hz and kept there for several minutes. Then the system is brought to rest. The 
cylinders are then brought to their initial rotational velocities (with the chosen 1/Ro), 
corresponding to a velocity for which the torque is accurate enough; generally of order 
2-3 Hz. We then slowly increase both velocities over 3-6 hours to their final velocities 
while maintaining 1/Ro fixed during the entire experiment. During this velocity ramp 
we continuously acquire the torque of this quasi-stationary state. The calibration of the 
system is done in a similar way; first we apply the maximum load on the system, going 
back to zero load, and then gradually adding weight while recording the torque. These 
procedures ensure that hysteresis effects are kept to a minimum, and that the system is 
always b rought to the s a me sta te before measuring. More details about the setup can be 
found in Ivan Gils et all (|2011al ). 

Local velocity measurements are done by laser Doppler anemomctry (LDA). We mea- 
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sure the azimuthal velocity component by focusing two beams in the radial-azimuthal 
plane. We correct for curvature effects of the outer cylinder by using a ray-tracer, see 
Huisman et all (|2012al ). The velocities are measured at midhcight (z = L/2) unless spec- 
ified otherwise. For every measurement position we measured long enough such as to 
have a statistically stationary result. 



4. Global response: Torque 

In this section, the global response of the TC system for the four simulated radius ratios 
is presented. This is done by measuring the scaling law(s) of the non-dimensional torque 
Nuu as function(s) of Ta. The transition between different types of local scaling laws 



in diff erent Ta-ranges is investigated, and rel ated to previous simulations (lOstilla et al. 
2013M ) and experiments ( van Gils et aZ.ll2012f) . 



4.1. Pure inner cylinder rotation 

The global response of the system is quantified by Nu u . By definition, for purely az- 
imuthal laminar flow, Nu u = 1. Once the flow is driven stronger than a certain critical 
Ta, large azimuthal roll structures appear, which enhance angular transport through a 
large scale wind. 

Figure [3] shows the response of the system for increasing Ta in case of pure inner 
cylinder rotation for four values of r\. Experimental and numerical results are shown in 
the same panels, convcring different ranges and t hus complementary, b ut consistent with 
each other. Numerical results for rj = 0.714 from lOstilla et all ( 2013a ) have been added 
to both panels. 



As has already been noticed in lOstilla et all (J2013M ) for r/ = 0.714, a change in the 
local scaling law relating Ta to Nu^ occurs at around Ta « 3 • 10 6 . We can interpret 



these changes in the same way as lOstilla et all (|20136[ ) and relate the transition in the 
Ta-Nuu, local scaling law to the break up of coherent structures. 

It is also worth mentioning that the exponent in the local scaling laws in the regime 
before the transition depends on the radius ratio. This can be seen in the compensated 
plot, and explains the curve crossings that we see in the graphs. 

For experiments (solid lines of figure [4}, a different local scaling law can be seen. 
In this case the experiments are performed at much higher Ta than the simulations. 
The scaling Nu^ ~ Ta 038 can be related to the so called " ultimate" regime, a regime 
wher e the boundary layers have become completely turbulent ( Grossmann fc Lohsell2011 > 



2012t Irluisman et all 120131) . As indicated for the case at n = 0.714 we expect that for 



increasing Ta also th e simulations become turbulent enough to reach this scaling law (cf. 
Ostilla et all (2013a))). In this regime, the local scaling law relating Ta and Nu^ has no 
dependence on n and thus is universal. 

In both experiment and simulation with the largest Ta, the value of r\ corresponding 
to the smallest gap, i.e. n = 0.909, has the highest angular velocity transport (Nu u ) at 
a given Ta. This can be phrased in terms of the pseudo-Prandtl-number a, introduced 
in EGL2007. As a smaller gap means a smaller a, we thus find a decrease of Nuu, 



with i ncreasing a, similarly as predicted iGrossmann 
(|2002l) for Nu(Pr) in RB convection for Pr > 1. 



Lohsej (|200lh and found Kia et 



4.2. Rossby number dependence 
In this subsection, the effect of outer cylinder rotation on angular ve locity transport will 



be studied. Previous experimental and numerical work at 77 = 0.714 (jPaoletti fc Lathrop 
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Figure 4: The global system response for pure inner cylinder rotation as function of the 
driving Ta: The top panel shows Nu u — 1 vs Ta for both simulations (points on the left 
of the graph ) and ex periments (lines on the right of the graph). Numerical data from 
Ostilla et all ( 2013al ) for 77 = 0.714 have been added to these figures. The bottom panel 
shows the compensated Nusselt ((Nu u — 1)/Ta 1 / 3 ) versus Ta. 



20 111 Ivan Gils et~al\ 12011 Zt lOstilla et al\ 1201361: iBrauckmann fc Eckhardl l2013al ) re- 
vealed the existence of an optimum transport where, for a given Ta, the transport of 
momentum is highest at a Rossby number Ro~ t , which depends on Ta and saturates 
around Ta ~ 10 10 . In this subsection, this work will be extended to the other values of 

V- 

Figure [5] shows the results of the numerical exploration of the Ro -1 parameter space 
between Ta = 4 ■ 10 4 and Ta = 2.5 • 10 7 . The shape of Nu^ = Nu UJ (Ro~ 1 ) curves and the 
position of Ro~p t depends very strongly on 77 in the Ta range studied in numerics. For the 
largest gap (i.e., 77 = 0.5), the optimum can be seen to be in the counter-rotating range 
(i.e., i?o _1 < 0) as long as Ta is high enough. On the other hand, for the smallest gap (i.e., 
77 = 0.909), the optimum is at co-rotation (i.e., Ro~ x > 0) in the whole region studied. 
The other values of studied 77 reveal an intermediate behavior. Optimum transport is 
located for co-rotation at lower values of Ta and slowly moves towards counter-rotation. 
For all values of 77, when the driving is increased, Ro~ t tends to shift to more negative 
values. 

The shift seen in the numerics may or may not continue with increasing Ta. The 
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Figure 5: Nu^ — 1 versus Ro^ 1 for the four values of rj studied numerically, 77 = 0.5 
(top), 0.714, 0.833 and 0.909 (bottom). The shape of the curve and the position of the 
maximum depend very strongly on both Ta and rj. 
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Figure 6: The left panel shows Nu^ versus Ta for all values of r\ and Ro~ l studied in 
experiments. The right panel shows the exponent 7 of the scaling law Nu^ ex Ta 1 for 
various Ro~ l , obtained by a least-square linear fit in log-log space. The average value of 
7 for each 77 is represented by the dashed lines, while the solid line represents the average 
value of 7 = 0.39 for all r\, which will be used for compensating Nu^. 



experiments conducted explore a parameter space of 10 10 < Ta < 10 13 and thus serve 
to explore the shift at higher driving. Figure [6] presents the obtained results. The left 
panel shows Nu u versus Ta for all measurements. The right panel shows the exponent 
7, obtained by fitting a least-square linear fit in the log-log plots. Across the n and 
Ro^ 1 range studied, the average exponent is 7 « 0.39. This value is used in figure [7] to 
compensate Nu u . The horizontality of all data points reflects the good scaling and the 
universality of this ultimate scaling behaviour Nu^ oc Ta 039 . 

Figure [8J shows the shape of the Nu^ = Nu ul (Ro^ 1 ) or alternatively Nu u = Nu UJ (a) 
curve at the cut-off region highlighted in figure [7] from experiment. The increased driving 
changes the characteristics of the flow. This is reflected in the very different shapes of 
the i?o _1 -dependence of Nu^ when comparing figures [5] and [51 and in the shift of Ro~p t . 

To summarize these effects, figure M presents both the 85% peak width ARo^ ax and 
the position of optimal transport Ro~^ t as a function of Ta and r\ obtained from numerics 



as well as the asymptotic value from experiments. The peak width Ai?o„ 
as: 



is defined 



ARo,„ nl , — 



Ro n 



Ro_ 



Nu^Ro-^dRo' 



(4.1) 



where Ro_q 85 and Ro g 5 



ma,x(Nu LU — 1) 

are the values of Ro~ x for which Nu M is 85% of the peak value. 
The 85% peak width can be seen to vary with driving, reflecting what is seen in figure 
[5] The shape of the Ro~ 1 -Nu u) curve is highly dependent of both r\ and Ta. Ro~^ t shows 
a very large variation across the Ta range studied in numerics. Therefore, one may ask 
the question: has the value of Ro~lt already saturated in our experiments? Figure [7] 
shows the trend for Nu^ for increasing Ta. This trend does not seem to vary much for 
different values of Ro~ x . Therefore, we expect the value of Ro^l to have already reached 
saturation in our experiments. 



We can c ompare these new experim e ntal res ults to the speculations made in lvan Gils et al 
( 20121 ) and Brauckmann fc EckhardtJ (J20136I ) for the dependence of the saturated Ro~^ t 
on n. This is shown in figure [101 Both speculated dependencies are shown to deviate 
substantially from the experimental results obtained in the present work. Even if the 
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Figure 8: The panels show Nu^/TaP' 39 versus either Ro~ l (top) or a (bottom) at the 
cut-off region highlighted in figure [7] for the values of rj studied experimentally. There is a 
strong 77-dependence of the curve Nu^/Ta ' 39 versus Ro~ x , even at the largest drivings 
studied in experiments. Optimal transport is located at Ro~l. = —0.20 for rj = 0.714, 
= -0.15 for r] = 0.769, Ro~l = -0.10 for 77 = 0.833 and Ro'l = -0.05 for 



Rt 



opt 



J opt 

0.33- 



77 = 0.909, corresponding toaw 0.33 — 0.35 for all values of 77. In the bottom panel, the 
maximum of the graph is less pronounced, i.e. it becomes more flat with increasing 77. In 
the limit rj — > 1, a does not tend to a finite limit, while Ro~ x does. This result highlights 
the advantage of using Ro~ x instead of o as a control parameter. 
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Figure 9: In the left panel, 85% peak width ARo^ ax vs Ta for the four values of 77 
analyzed in numerics. The peak width can be seen to vary with driving, and for smaller 
gaps is larger for larger values of Ta. In the right panel, Ro~ t vs Ta for the same four 
values of 77. The location of the optimal transport has a very strong dependence on the 
driving, especially for the largest values of rj. As driving increases beyond the numerically 
studied range and overlaps with experiments, Ro~ pt should tend to the experimentally 
found values, represented as dashed lines in the figure. The asymptote for 77 = 0.5 is 
obtained from iMerbold et all (|2013l ). The trend appears to be less clear for i] = 0.909, 
but this might be understandable from the peak width at the highest driving Ta. 
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speculation from Ivan Gils et al\ (J2012J) appears to be better for this 77-range, for previous 
experimental data at rj — 0.5, it is in clear conflict with the experimentally measured 
value for optimal transport by IMerbold et all ( 2013 ) . 

This section has shown that the radius ratio has a very strong effect on the global 
response and especially on optimal transport. Significantly increased transport for co- 
rotation has been found at the lowes t drivin gs based on the DNS results. This finding 
was already reported in lOstilla et all (J20136I ) for ?? = 0.714, but the transport increase 
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was marginal. For n = 0.833 and especially for n = 0.909 the transport can be increased 
up to three times. The shift of Ro~ t has also been seen to be much bigger and to happen 
in a much slower way for smaller gaps. The reason for this will be studied in Section 
using the local data obtained from experiments and numerics. 



5. Local results 

In this section, the local angular velocity profiles will be analyzed. Angular velocity is 
the transported quantity in TC flow and shows the interplay between the bulk, where 
the transport is convection dominated, and the boundary-layers, where the transport is 
diffusion dominated. Numerical profiles and experimental profiles obtained from LDA 
will be shown. The angular velocity gradient in the bulk will be analyzed and connected 
to the optiiual transport. In addition, the boundary layers will be analyzed and compared 
to the results from the analytical formula from EGL 2007 for the BL thickness ratio in 
the non-ultimate regime. 

5.1. Angular velocity profiles 

Angular velocity u> profiles obtained from numerics are shown in figure 1111 Results are 
presented for four values of r\ and selected values of Ro^ 1 at Ta = 2.5 • 10 7 (and Ta = 
2.39 • 10 7 for r\ = 0.714). Experimental data obtained by using LDA are shown in figure 
[l2]for three values of r\: from top-left to bottom, ij = 0.714 for Rei — Re = 10 6 , rj = 0.833 
for Ta = 5 • 10 11 , and r\ = 0.909 for Ta = 1.1 • 10 11 . 

The different radius ratios affect the angular velocity profiles on both boundary layers, 
as the two boundary layers are more asymmetric for the wide gaps; and they affect the 
bulk, as the bulk angular velocity is smaller for wide gaps. These effects will be analyzed 
in the next sections. 

5.2. Angular velocity profiles in the bulk 

We now analyze the properties of the angular velocity profiles in the bulk. We find that 
the slope of the profiles in the bulk is controlled mainly by Ro~ x and less so by Ta. This 
can be understood as follows: The Taylor number Ta acts through the viscous term, 
dominant in the boundary layers, while Ro~ l acts throug h the Coriolis force, p resent in 



the whole domain. These results extend the finding from lOstilla et al\ (J2013M ) to other 
values of n. 

To further quantify the effect of Ro^ 1 on the bulk profiles, we calculate the gradient 
of (u)) z . This is done by numerically fitting a tangent line to the profile at the inflection 
point using the two neighboring points on both sides; such fit is shown in the left panel 
of figure [j"3l 

Figure [14] shows four panels, each containing the angular velocity gradient in the bulk 
from the numerical simulations and experiments for a given value of r\. We first notice 
that the angular velocity gradients from experiment and numerics are in excellent agree- 
ment. Next the connection between a fiat angular velocity profile and optimal transport 
can n ow be seen for other values of r\ and not just for r\ = 0.714 as reported previ- 
ously (Ivan Gils et alJl2011a ). Once Ro~ x < Ro~p t , the large scale balance analyzed in 



Ostilla et all (J2013M ) breaks down, and a "neutral" surface which reduces the transport 



appears in the flow. 

In simulations, because of resolution requireme nts, we are unable of driving the flow 



strongly enough to see a totally flat bulk profile. In lOstilla et all (J20136I ). a linear extrap- 
olation of the bulk angular velocity gradient was done to give an estimate for the case 
when this profile would become horizontal, i.e., d{uj) z /df = 0, and thus give an estimate 
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Figure 11: Aziniuthally, axially and temporally averaged angular velocity (w) z versus 
radius r for: r) = 0.5, r) = 0.714, t] = 0.833, and 77 = 0.909. Data is for Ta = 2.5 ■ 10 7 
(Ta = 2.39 • 10 7 for 77 = 0.714) and selected values of Ro~ x . For smaller 77, the w-bulk 
profiles differ more from a straight line, and have, on average, a smaller value. 



of Ro opt . For 77 = 0.714 this estimate agreed with the numerical result within error bars. 



Here, w e extend this analysi s for the other values of 77 and, as we shall see, successfully. 
Ostilla et al\ ( 2013a ). an almost linear relationship between Ro~ x and d{w) z /df 
seen. This linear relationship is extrapolated and plotted in each panel. This 



can 



As m 
be 



extrapolation gives an estimate for Ro opt (Ta — > 00), which we can compare to the ex- 
perimentally determined Ro~ p \(Ta —> 00). For ?] = 0.833, Ro~ 1 (Ta — > 00) s=s —0.12 
corresponding to a w 0.38 is obtained, and for 77 = 0.909, Ro~ 1 (Ta — > 00) w —0.05, 
corresponding to a « 0.31 is obtained. These values are (within error bars) also obtained 
for Ro~ t at the large Ta investigated in experiments, namely Ro~ t = —0.10 and —0.05, 



respectively. 

For 77 = 0.5, Ro~ pt (Ta — > 00) « —0.33 is obtained, corresponding to a » 0.2. This value 
is consistent with the numerical results in lBrauckmann fc Eckhardtl(|2013M ). which report 
o-opt ~ 0.2. However, care must be taken, as fitting straight lines to the w-profiles gives 
higher residuals for 77 = 0.5 as the profiles deviate from straight lines (cf. top-left panel 
of figure ITT|) . A fit to the "quarter-Couette" profile derived from upper bound theory 
( Busselll967l ) is much more approp riate for 77 = 0.5 at the strongest drivings achieved in 
experiments (JMerbold et aZ.ll2013i ). This is because the flow feels much more the effect 
of the curvature at the small 77. On the other end of the scale, the linear relationship 
works best for smallest gaps, i.e. r\ = 0.909 (cf. the bottom right panel of figure fTTj) where 
curvature plays a small effect. 



Optimal Taylor- Couette flow: Radius ratio dependence 



19 




-•- 


1/Ro=0.34 


-*- 


0.20 


-♦- 


0.09 


— •— 


0.00 


A 


-0.07 


-«- 


-0.13 


-*- 


-0.18 


-k- 


-0.20 


-m- 


-0.23 


-»- 


-0.27 


-♦- 


-0.30 


-*— 


-0.33 


-A- 


-0.40 


-T- 


-0.53 




Figure 12: Angular velocity profiles obtained by LDA for n = 0.714 at either Rei — 
Re = 10 6 (top-left), n = 0.833 or at Ta = 5 • 10 11 (top-right) and r\ = 0.909 at Ta = 
1.1- 10 11 (bottom), to explore different dependencies in parameter space. Data is taken 
at a fixed axial height, but as the Taylor number Ta is much larger than in the numerics, 
the axial dependence is much weaker. 



5.3. Angular velocity profiles in the boundary layers in the classical turbulent regime 

As the driving is increased, the transport is enhanced. To accommodate for this, the 
boundary layers (BLs) become thinner and therefore the w-slopes (d r uj) become steeper. 
Due to the geometry of the TC system an intrinsic asymmetry in the BL layer widths is 
present. More precisely, the exact relationship d r (cj)\ = ?fd r (uj)\i holds for the slopes 
of the boundary layers, due to the r independence of J w , cf. EGL 2007 and eq. (|2.5l) . 

An analysis of the boundary layers was not possible in the present experiments because 
the present LDA measurements have insufficient spatial resolution to resolve the flow in 
the near-wall region. Therefore, only DNS results will be analyzed here. In simulations 
the driving is not as large as in experiments, and as a consequence the shear in the BLs 
is expected to not be large enough to cause a shear-instability. This means that the BLs 
are expected to be of Prandtl-Blasius (i.e. laminar) type, even if the bulk is turbulent. 
On the other hand, in the experiments both boundary layers and bulk arc turbulent, i.e. 
the system is in the "ultimate regime" . 

Using the DNS data, we can compare the ratio of the numerically obtained boundary 
layer widths with the analytical formula for this ratio obtained by EGL 2007 for laminar- 
boundary layers, namely: 






\W — UJbulk] 



(5.1) 
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Figure 13: An example of the two fitting procedures for the bulk angular velocity gradi- 
ent and for boundary layer thicknesses is shown here. Both panels show the 9,z, and t 
averaged azimuthal velocity and angular velocity for r\ = 0.5, Ta=\ ■ 10 7 , and pure inner 
cylinder rotation. In the left panel a line is fitted to the bulk of the angular velocity to 
obtain the bulk gradient. The right panel shows the three- lines-fit to the whole profile 
to obtain the width of its boundary layers, used in Section I5~3l Both bulk fits are done 
at the inflection point, but for different variables (a) or ug), which gives slightly different 
slopes (and intersection points). 



where the value of Wbulk is some appropriate value in between for which the angular 
velocity at the inflection point of the profile might be chosen, i.e., the point at which the 
linear bulk profile fit was done to obtain A° and A^. The value uibuik is taken from the 
numerics, and may bias the estimate. 

To calculate the boundary layer thicknesses, the profile of the mean azimuthal velocity 
(ug) z is approximated by three straight lines, one for each boundary layer and one for the 
bulk. For the boundary layers the slope of the fit is calculated by fitting (by least-mean- 
squares) a line through the first two computational grid points. For the bulk, first the 
line is forced to pass through the grid point which is numerically closest to the inflection 
point of the profile. Then its slope is taken from a least mean square fit using two grid 
points on both sides of this inflection point. The respective boundary layer line will cross 
with this bulk line at a point which then defines the thickness of that boundary layer. 

The results obtained for A°/A^ both from equation 15. II and directly from the simula- 
tions is shown in figure Q3J Results are presented for the four values of r\ and only for 
the highest value of Ta achieved in the simulations. The boundary layer asymmetry for 
counter-rotating cylinders (i.e., Ro" 1 < 0) grows with larger gaps. This is to be expected, 
as the rj~ 3 term is much larger (s=s 8) for the largest gap as compared to the smallest gap 
(~ 1.3). This is consistent with the ij and thus a restriction in EGL 2007 to a range of 
smaller gap widths. 

As noticed already in lOstilla et al\ (20136) we find that the fit is not satisfactory for 
co-rotation (i.e., Ro^ 1 > 0) at the lowest values of r], but is satisfactory for counter- 
rotation (ie. Ro" 1 < 0). This can be explained from the fact that in EGL 2007, equation 
(|5.1[) is obtained by approximating the profile by three straight lines, two for the BLs 
and a constant w line for the bulk. This does not hold very well for co-rotating cylinders, 
where the bulk has a steep gradient (see figure IT~5]) . The steepness of the gradient in the 
bulk is larger for smaller values of r\ as seen in Section 15.21 This is why the formula fails 
for r\ = 0.5. 
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Figure 14: Bulk angular velocity gradient d(uj) z /df against Ro -1 for the four values of rj 
explored in simulations, rj = 0.5 (top left), 77 = 0.714, (top right), 77 = 0.833, (bottom left), 
and 77 = 0.909 (bottom right). Data from experiments obtained by LDA is also plotted 
for the three values of 77 for which it was experimentally measured (green circles). For 
all values of 77 except 77 = 0.5, for co-rotation and slight counter-rotation there is once 
again an almost linear relationship between Ro~ l and d(uj) z /dr. A black straight line is 
added to extrapolate this relationship in order to estimate Ro~ t . A plateau, in which 
the radial gradient of (ui) z is small can be seen around optimal transport, indicating a 
large convective transport of angular velocity 

For co-rotation the boundary layers are approximately of the same size, and the ratio 
A° /AJ, is very close to 1. If one inverts equation (|5.ip by aproximating this ratio by 1, 
an estimate of what the angular velocity will be in the bulk due to the boundary layer 
slope asymmetry is obtained: 



Ubulk 



corresponding to: 



V 



OJ 



bidk 



1 



3 P 

rfw\ 



V 



(5.2) 



(5.3) 



in the lab frame. This expression gives an estimate for Ubuik when the profile is flattest, 
and has been represented graphically in figure [TT] Indeed, one can take this estimate (for 
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Figure 15: A°/A^ from simulations (dots) and from equation (|5.1[) (solid lines) versus 
Ro^ 1 for the four values of n studied numerically, i] = 0.5 (top-left), 0.714 (top-right), 
0.833 (bottom-left) and 0.909 (bottom-right) at Ta = 2.5 • 10 7 . The numerical results and 
the estimate from equation (|5.ip match very well for larger values of r] and especially for 
counter-rotating cylinders (l/Ro < 0). The asymmetry between the boundary layers can 
be seen to be larger for smaller values of n, which is expected as equation (|5.1|) contains 
the explicit factor n~ 3 . 



example 0.27 for n = 0.714) and compare it with figures ITT1 and [j~2l We note that the 
value of uj in the bulk for the flattest profile in the numerics (at Ro^ 1 w Ro^ t {Ta)) lies 
around uJbuik- We can also notice that the profiles for Ro" 1 > Ro~p t approximately cross 
each other at the same point, and this point has a value of ui ~ uibuik- This effect can only 
be seen in the numerics, as these approximations break down once the boundary layers 
become turbulent. The cross points of the curves are taken as an estimate for LObuik-, and 
this is represented against eq. 15.21 in figure [TBI 

To understand why the boundary layers are of approximately the same thickness de- 
spite the different initial slopes at the cylinders one has to go back to equation (|2.5|) . The 
angular velocity current has a diffusive part and a convective part. Per definition in the 
boundary layer the diffusion dominates and in the bulk the convection does. Thus the 
boundary layer ceases when convection becomes significant. But convection is controlled 
by the wind. Thus in essence the boundary layer size is controlled by the wind and not 
immediately by the initial slope at the wall. Due to continuity, if the rolls penetrate the 
whole domain the wind may be expected to be the same close to the inner and close to 
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Figure 16: uJbuik as a function of 77, taken from both eq. 15.21 and from the crosspoints of 
the w-profiles in figure QTJ The trend is the same in both data sets. A smaller value of r\ 
decreases the value of the bulk angular velocity. 



the outer cylinder. This suggests that the flow organizes itself in a way that the boundary 
layer extensions (or widths) might be similar, even if the initial slopes at the walls are 
different. 

What happens for counter-rotation, or more precisely when Ro~ l < Ro~ t l For Ro~ x 
below the optimum Ro~ t a so-called neutral surface will be present in the flow, which 
separates the Rayleig h-stable and -unstab le areas. The wind drastically changes in the 
Rayleigh-stable areas (|Ostilla et al\\2013m , leading to very different wind velocities close 
to the outer and inner cylinder, respectively. The wind at the outer cylinder will be 
weaker, as the rolls cannot fully penetrate the Rayleigh-stable domain. This means that 
the outer cylinder boundary layer will extend deeper into the flow, in accordance to what 
is seen in figure [T51 



6. Summary and conclusions 

Experiments and direct numerical simulations (DNS) were analysed to explore the 
effects of the radius ratio r\ on turbulent Taylor-Coucttc flow. Numerical results corre- 
sponding to Taylor numbers in the range of 10 4 < Ta < 10 8 alongside with experiments 
in a Taylor number range of 10 10 < Ta < 10 12 were presented for four values of the 
radius ratio rj. 

First the influence of the radius ratio on the global scaling laws Nu u ~ Ta 1 was 
studied. The local scaling exponent 7 describing the response of the torque caused by 
a Taylor number increase, is barely modified by varying the radius ratio ?y. Indeed, in 
experiments a universal exponent of 7 rj 0.39 is obtained, independent of radius ratio 
and outer cylinder rotation. For the numerical simulations at lower Ta similar universal 
behavior can be observed. The transition associated to the vanishing of coherent struc- 
tures can also be appreciated at Ta ~ 10 6 for all values of r\. Before this transition local 
exponents of 7 ~ 0.33 are seen and after the transition these decrease to about 7 « 0.2. 

The radius ratio does play a very important role in optimal transport. At smaller gaps, 
i.e., for larger rj, at the lower end of the Ta range a very large increase in transport for 
corotating cylinders can be seen. The shift towards the asymptotic optimal transport 
happens in a much slower way for small gaps, but this shift is seen for all studied radius 
ratios. For the largest gap (rj = 0.5), optimal transport for pure inner cylinder rotation at 
the lowest drivings is obtained. The shift towards the asymptotic value happens suddenly, 
as two peaks can be seen in the Nu^ versus Ro~ x curve, and one of the peaks becomes 
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larger than the other one as driving increases. This might point in the direction of 
different phenomena happening at larger gaps. Finally, the a symptotic values of R o~ f 



obtained in experiments were com pared to the speculations of lvan Gils et al\ (J2012I ) and 



Brauckmann fc Eckhardtl (J2013M ). Both of these speculations were found to deviate from 



experimental and numerical results. 



When looking at the local results, as in lOstilla et al\ (12013a) we can link the optimal 



transport in the smallest gaps to a balance between Coriolis forces and the inertia terms 
in the equations of motion. The flattest profiles in the bulk are linke d to optimal transport 



in experiments. With the numerics the extrapolation presented in lOstilla et all (|20136[ ) 
for predicting optimal transport was extended to other radius ratios. It is found to work 
well for all selected r\ except for rj = 0.5. At this r), i.e., for the largest gap considered 
here, the most obvious problem is that the ui profiles strongly feel the effect of curvature 
difference at the inner and outer cylinders and a straight line fit to the bulk is not 
appropriate. There may be additional reasons for this discrepancy and optimal transport 
in large gaps requires more investigation. 

Finally, the boundary layers have been analyzed. The outer boundary layer is found 
to be much thicker than the inner boundary layer when Ro^ 1 < Ro~ t . We attribute 
this to the appearance of Rayleigh- stable zones in the flow. This prevents the rolls from 
covering the full domain between the cylinders. As the boundary layer size is essentially 
determined by the wind, if the rolls penetrate the whole domain (which is the case for 
Ro^ 1 > Ro~ t ), both boundary layers arc approximately of the same size. If the rolls do 
not penetrate the full domain, the outer boundary layer will be much larger than the 
inner boundary layer, in accordance with the smaller initial slope of ui(r) at the cylinder 
walls. 

In this work, simulations and experiments have been performed on a range of radius 
ratios between 0.5 ^ n ^ 0.909. Insi g hts for the small gaps seem to be consistent with 



what was discussed in lOstilla et al\ (J2013M ). However, for rj = 0.5 the phenomena of 



optimal transport appears to be quite different. Therefore, our ambition is to extend the 
DNS towards values of rj smaller than 0.5 to improve the understanding of that regime. 
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